
#include "population.h"
#include <math.h>
#include <iostream.h>
#include <iomanip.h>

population::population() {
	popSize = 0;
}

population::~population( ) {}

population::population( const population& rhs ) {
	popSize = rhs.popSize;
	pop.resize(popSize);

	for(int i=0; i<popSize; i++) {
		pop[i] = rhs.pop[i];
	}
}

population::population( int aSize  ) {
	popSize = aSize;
	pop.resize(popSize);
	
}

population::population( string& popFileName, string& molSettingFileName,  vector<double>(*anEvaluator)(vector<atom>& a, int optimize) ) {
	
	ifstream popFile( popFileName.c_str() );
	
	string dummy;
	int nAtoms;
	popFile >> dummy >> popSize >> dummy >> nAtoms;

	pop.resize(popSize);

	for(int i=0; i<popSize; i++) {
		pop[i] = molCluster( molSettingFileName, anEvaluator ); 
		for( int j=0; j<nAtoms; j++ ) {
			popFile >> pop[i].atoms[j].label >> pop[i].atoms[j].x >> pop[i].atoms[j].y >> pop[i].atoms[j].z;
		}
	}

	popFile.close();


}


void population::selectOneBest( molCluster& aMolCluster ) {
	sort();
	aMolCluster = pop[0];
}

void population::selectNBest( int n, vector<molCluster>& molClusters ) {
	sort();
	for(int i=0; i<n; i++) {
		molClusters[i] = pop[i];
	}
}
        

void population::selectOneRandom( molCluster& aMolCluster ) {
	
	aMolCluster = pop[ (int)( (int)( mt()*10000 )%popSize) ];

}

void population::selectNRandom( int n, vector<molCluster>& molClusters ) {
	
	for(int i=0; i<n; i++) {
		selectOneRandom( molClusters[i] );	
	}

}

void population::selectOneWorst( molCluster& aMolCluster ) {
	sort();
	aMolCluster = pop[popSize-1];
}

void population::selectNWorst( int n, vector<molCluster>& molClusters ) {
	sort();
	for(int i=0; i<n; i++) {
		molClusters[i] = pop[popSize-n+i];
	}
}

        
void population::selectOneTournament( int tournamentSize, molCluster& aMolCluster ) {
	
	vector<molCluster> molClusters(tournamentSize);
	for(int i=0; i<tournamentSize; i++) selectOneRandom( molClusters[i] );

	aMolCluster = molClusters[0];
	int bestIndex = 0;
	for(int i=1; i<tournamentSize; i++) {
		if( molClusters[i].properties[0]< molClusters[bestIndex].properties[0] ) {
			bestIndex = i;
			aMolCluster = molClusters[i];
		}

	}

}

void population::selectNTournament( int n, int tournamentSize, vector<molCluster>& molClusters ) {

	for(int i=0; i<n; i++) {
		selectOneTournament( tournamentSize, molClusters[i]);
	}	

}
        
void population::selectNSUS( int n, vector<molCluster>& molClusters ) {

	sort();
	int step = popSize/n;

	int start = (int)(( mt()*10000 ))%popSize;
	molClusters[0] = pop[start];
	int nSelected = 1;

	while(nSelected<n) {
		molClusters[ nSelected ] = pop[ (start+step*nSelected)%popSize ];	
		nSelected++;
	}

}

void population::replaceMuLambda( int nElitist, population& rhs ) {

	sort();
	rhs.sort();

	for( int i=nElitist; i<popSize; i++ ) {
		pop[i] = rhs[ i-nElitist ];
	}

}

molCluster& population::operator[]( int index ) {
	
	return pop[index];
	
}
        
population& population::operator=( const population& rhs ) {

	popSize = rhs.popSize;
	pop.resize(popSize);
	for(int i=0; i<popSize; i++) {
		pop[i] = rhs.pop[i];
	}
	return (*this);
}

population& population::operator+( const population& rhs ) {

	population newPop( rhs.popSize + popSize );
	int i;

	for(i=0; i<popSize; i++) {
		newPop[i] = pop[i];
	}


	for(i=popSize; i<rhs.popSize+popSize; i++) {
		newPop[i] = rhs.pop[i-popSize];
	}

	return (*this) = newPop;

}
        

population& population::operator+=( const population& rhs ) {

	return (*this)+rhs;

}

void population::sort() {
	sortMolClusterVector( pop );
}
        
double population::bestFitness() {

        vector<double> array(popSize);

        int i;

        for(i=0; i<popSize; i++) {
                array[i] = pop[i].properties[0];
        }

        sortDoubleVector( array );

        double best;
        best = array[0];

        return best;

}
        
double population::worstFitness() {

        vector<double> array(popSize);

        int i;

        for(i=0; i<popSize; i++) {
                array[i] = pop[i].properties[0];
        }

        sortDoubleVector( array );

        double worst;
 	worst = array[popSize-1];

	return worst;
}

double population::medianFitness() {

        vector<double> array(popSize);

	int i;

        for(i=0; i<popSize; i++) {
                array[i] = pop[i].properties[0];
        }

	sortDoubleVector( array );

	double median;
	if(popSize%2) median = array[ (int)floor(popSize/2) ];
	else median = 0.5 * ( array[ (popSize/2) ] + array[ (popSize/2)-1 ] );

	return median;
}
        
double population::meanFitness() {
	vector<double> array(popSize);

        double average = 0;
        for(int i=0; i<popSize; i++) {
                array[i] = pop[i].properties[0];
                average += array[i];
        }
        return average/=popSize;

}

double population::stdevFitness() {

	vector<double> array(popSize);
	int i;

	double average = 0;
	for(i=0; i<popSize; i++) {
		array[i] = pop[i].properties[0];
		average += array[i];
	}
	average/=popSize;

	double sum=0;
	for(i=0; i<popSize; i++) {
		sum += (array[i]-average)*(array[i]-average);
	}
	
	return sqrt( sum/(popSize-1) );

}

int population::size() {
	return popSize;
}
	
void population::printToFileStream( ofstream& out ) {

        //output size
        out << scientific << setprecision( PRECISION );
	out << "popSize	" << size() << endl;
	out << "nAtoms	" << pop[0].size() << endl << endl;

	for(int i=0; i<popSize; i++) out << pop[i];
}


ostream& operator<<( ostream& out, population& aPop ) {
        
	//output size
        out << scientific << setprecision( PRECISION );
        out << "popSize " << aPop.size() << endl;
        out << "nAtoms  " << aPop.pop[0].size() << endl << endl;

        for(int i=0; i<aPop.size(); i++) out << aPop.pop[i];

	return out;
}

